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Abstract. 

For optimal processing and design of entangled polymeric materials it is important 
to establish a rigorous link between the detailed molecular composition of the polymer 
and the viscoelastic properties of the macroscopic melt. We review current and past 
computer simulation techniques and critically assess their ability to provide such a link 
between chemistry and rheology. We distinguish between two classes of coarse-graining 
levels, which we term coarse-grained molecular dynamics (CGMD) and coarse-grained 
stochastic dynamics (CGSD). In CGMD the coarse-grained beads are still relatively 
hard, thus automatically preventing bond crossing. This also implies an upper limit on 
the number of atoms that can be lumped together (up to five backbone carbon atoms) 
and therefore on the longest chain lengths that can be studied. To reach a higher 
degree of coarse-graining, in CGSD many more atoms are lumped together (more than 
ten backbone carbon atoms), leading to relatively soft beads. In that case friction 
and stochastic forces dominate the interactions, and actions must be undertaken to 
prevent bond crossing. We also review alternative methods that make use of the tube 
model of polymer dynamics, by obtaining the entanglement characteristics through a 
primitive path analysis and by simulation of a primitive chain network. We finally 
review super-coarse-grained methods in which an entire polymer is represented by a 
single particle, and comment on ways to include memory effects and transient forces. 

PACS numbers: 66.20.Cy,47.57.Ng,83.80.Sg 
1. Introduction 

Ever since tlie first synthetic polymers were made by Staudinger and their molecular 
weights were measured by Debye and Bueche, polymer systems have received continual 
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interest from both theorists and experimentahsts. Most processing of polymers takes 
place in the melt state, when they are very viscous and have surprising (temporary) 
elastic properties. For optimal processing, design, and application of polymeric melts 
it is important to establish a rigorous link between their molecular composition and 
macroscopic mechanical properties, i.e. between polymer chemistry and rheology. 

Establishing such a link is by no means an easy task and has despite considerable 
progress not yet been fully achieved. Unlike simple molecules, interactions within 
and between polymer chains are characterised by a large range of different time and 
length scales. The dynamics of atomic bond vibrations are characterised by Angstrom 
and sub-picosecond scales, whereas the dynamics of statistically independent Kuhn 
segments are characterised by nanometer and tens of picosecond scales [I] . Beyond this, 
the connectivity and mutual uncrossability of the segments cause an interdependence 
between features on a hierarchy of scales. The full chain, finally, is characterised by its 
radius of gyration of 10 to 100 nanometers. This may not appear as a very large increase 
from the Kuhn segment scale, but the associated time scales increase dramatically, to 
milliseconds, seconds or even longer. 

The enormously long intrinsic time scales are often rationalized by viewing polymer 
systems as temporary rubbery networks. Such a network arises as a result of mutual 
uncrossability of the polymer chains - they are entangled. With the advent of reptation 
theory of de Gennes, Doi and Edwards |1;|, a new concept was introduced in the theory 
of polymer dynamics. In reptation theory each polymer is supposed to move in a tube 
around a Gaussian path in space. The tube only serves one purpose, namely to roughly 
represent the uncrossability of the surrounding chains and to turn a difficult multichain 
problem into a one-chain problem. The tube clearly is an idealised mean field concept, 
and many extensions of the original model have been postulated to explain experimental 
observations; the difficulty is that without detailed information such postulates often 
remain unchecked [2]. 

In this review we focus on computer simulation techniques that may help 
establishing the link between chemistry and large scale dynamics and rheology. 
Computer simulations can check the assumptions made in theoretical models, and 
either accept or reject them [2]. Ideally, one would like to use atomistically detailed 
models because they can accurately capture the subtleties in the interactions between 
the polymers. However, because of the intrinsically slow dynamics, a direct prediction 
of the large scale dynamics and rheology from such detailed models is computationally 
unattainable (except for the case of unentangled or slightly entangled polymers). 
Traversing the road from chemistry to rheology is usually a bumpy ride: the only way to 
reach the larger scales is by going through a succession of coarse-graining steps. Different 
levels of coarse-graining are schematically indicated in Fig. [T| 

A number of excellent reviews on coarse-graining of entangled polymers have 
already appeared in the literature. Some of these reviews focus mainly on static 
properties [^-[6], whereas others focus on the dynamic properties as well j^js]. This 
review provides an up-to-date overview of past and current coarse-graining efforts, with 
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Figure 1. The dynamics of polymers can be simulated by different methods which 
differ in their level of detail describing polymer configurations and entanglement effects. 
Top: atomistic molecular dynamics simulations are usually sufficiently detailed to 
faithfully predict the dynamic properties of polymers of a specified chemistry. Second 
from top: small groups of atoms are lumped together (up to five backbone carbon 
atoms) into coarse-grained entities; the interaction between these 'beads' are still 
hard enough to make bond crossing energetically unfavourable. Middle: many more 
atoms are lumped together (more than 10 backbone carbon atoms) into coarse-grained 
entities; the interactions between these 'blobs' are so soft that uncrossability needs to 
be enforced by additional constraints. Second from bottom: polymer chains are forced 
to diffuse along primitive paths, representing the temporary topological network which 
arises because of uncrossability of chains. Bottom: each polymer chain is represented 
by just one particle; the entanglement effect is captured through addition of slow 
variables whose deviations from equilibrium lead to transient forces. 



a strong emphasis on our own opinion about the correct way to proceed in the future. 
We therefore focus on methods to handle friction, uncrossabihty and memory, all of 
which are intimately connected to the long time dynamic properties and rheology. 

The review is structured as follows. In section |2] we review the work that has 
been done using atomistically detailed models. In section |3] we introduce the theory 
of coarse-graining and make a classification of coarse-grained dynamics methods into 
two classes: 1) those with relatively low levels of coarse-graining (up to five backbone 
carbon atoms), where friction is not dominant or ignored and where uncrossability is 
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included automatically by the excluded volume interactions between relatively hard 
beads, and 2) those with relatively high levels of coarse-graining (more than ten 
backbone carbon atoms), where friction is dominant and measures need to be taken 
to ensure uncrossability. Section |4] reviews the work that has been undertaken in the 
first class and section |5] reviews the work done in the second class. In section [6] we 
review the attempts that have been undertaken at characterising primitive paths and 
the use of primitive primitive path models. In section [7] we focus on the ultimate form of 
particulate coarse-graining, where an entire polymer is represented by just one particle. 
Finally, we conclude in section [8j 

2. Atomistically detailed models 

2.1. Possibilities and limitations 

On a microscopic level, molecular dynamics (MD) and Monte Carlo (MC) simulations 
can be performed, in which each atom of a polymer chain is represented separately, see 
Fig.|l|(top) (To). The atoms are modelled as interacting particles and their positions are 
updated according to Newton's laws (MD) or by certain trial moves which are accepted, 
depending on the statistical ensemble, with a certain probability (MC). Accurate force 
fields have been constructed to cater for anyone's research interests, provided they 
exclude chemical reactions and other phenomena of a quantum mechanical nature. Bulk 
behaviour is simulated by applying periodic boundary conditions to the simulation box. 
Typical MD simulations cover the motion of a few tens of thousands of atoms over a 
period of a few nanoseconds; on current computers such a run would take a week to 
complete. 

This already sets a limit to the length of polymers that can be simulated in atomistic 
detail. First, the radius of gyration of the polymer chain typically grows as ^/n., where 
n is the number of monomers per chain. A polymer should not interact with itself 
via the periodic boundaries, which means that the volume of the box, and hence the 
number of particles, should scale as n^/^. Second, the longest intrinsic relaxation time 
of a polymer chain in a melt scales very fast, usually as for unentangled polymers 
and as n^'^ for entangled polymers [l]. To measure long-time correlation functions, the 
simulations must be performed at least as long as this time scale. 

Sometimes one is interested in the static properties of a melt of long polymer 
chains, or in faster dynamical relaxation processes on a more local scale. In such cases, 
one still needs to make sure that the atomistic conformations in the simulation box have 
sufficiently relaxed to faithfully represent an equilibrium state. 

2.2. Equilibration schemes for atomistically detailed models 

So how does one ensure that an atomistically detailed polymer melt is sufficiently 
relaxed? A first reaction may be to start with an ensemble of chains with a correct 
end-to-end distance distribution, arrange them randomly in the simulation box, and 
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Figure 2. (color online) Schematic of the end-bridging (EB) move. Left: initial 
configuration. The red arrow indicates that atom i of chain ich will attack atom j 
of chain jch] then, trimer {ja,jb,jc) will be excised from the system. Right: chain 
configurations after the application of the move. Atoms i and j are connected through 



the new trimer bridge (ja' , jb' , jc'). Reproduced from 19 



introduce excluded volume rapidly. However, Auhl et al. 11 showed that this procedure 
leads to deformations on short length scales, which relax only when the chains have 
moved over their own size, i.e. after one longest intrinsic relaxation time. The authors 
also showed how this local deformation may be overcome by first pre-packing the 
Gaussian chains, reducing the density fluctuations in the system, followed by a more 
gradual introduction of excluded volume. Another method is to apply a double-bridging 
Monte Carlo algorithm in which new bonds are formed across a pair of chains, creating 
two new chains substantially different from the original two [TT] . 

The latter method is similar in spirit to the end bridging (EB) Monte Carlo method, 
see Fig. |2| which has been developed extensively in the groups of Mavrantzas and 
Theodorou [l2 22 . A thorough analysis of the geometric formulation and numerical 
implementation of the original EB algorithm is given in Ref. 12 . Using this method. 



equilibrated samples of polyethylene (up to C500) [12] and cis-l,4-polyisoprene (up to 
C200) [14, 15 have been generated. Later faster methods were developed with directed 



moves, dubbed Directed Internal Bridging (DIB) and Directed End Bridging (DEB), 
which enabled them to equihbrate polyethylene samples as long as Ceooo [l6||l7] . 
Extensions of the end-bridging algorithm to different chain architectures have also been 



introduced 18 



An exciting development, of relevance to the rheology of polymer melts, is the 
use of end bridging moves in combination with a small tensorial field to obtain well- 
equilibrated, pre-oriented, strained configurations 13,20,21. This way, viscoelastic 



properties can be measured from atomistically detailed models, although up to now 
only simulations with unentangled chains have been reported. 

In light of section |4.3| we already mention that atomistic conformations can also 
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be equilibrated efficiently by mapping them onto a coarse-grained model, equilibrating 
the coarse-grained model, and finally reinserting the atomistic details. Of course End 



Bridging moves can also be applied in these coarse-grained equilibration runs 23 25 



2. 3. Examples of atomistic simulations of dynamic properties 

Given the very steep increase of terminal relaxation times with molecular weight, it 
should come as no surprise that atomistically detailed MD simulations of the dynamics 
of polymer melts have only been performed for relatively short unentangled or slightly 
entangled chains. 

Polyethylene, because of its simplicity, has been one of the most favourite polymers 
studied in the past. Starting in 1995, Paul and co-workers studied the diffusion 
coefficient of unentangled chains, C44 26 and Cioo 27 , both by an explicit atom 



model and by a united atom model. In this united atom model, which was optimised 
to reproduce experimental P-V-T behaviour, each methyl CH^; group is treated as one 
particle. The viscosity of Cioo was investigated by Moore et al. 28 . Later, Harmandaris 
et al. studied the diffusion coefficient of samples of C24, Cyg and Cisg polyethylene which 
had been previously relaxed by the EB algorithm 29 . The shear relaxation modulus 
G{t) of the C24 and Cyg samples were also obtained by analysing equilibrium stress 



fluctuations 13 . Soon after, in 2001, our group reported the mean square displacement, 
shear relaxation modulus, and viscosity of a C120 polyethylene sample 30 . Finally, in 



2003 mean square displacements and diffusion coefficients for chain lengths C78 to C 



250 



were reported 31 



Other polymer species have been studied as well, notably 1,4-polybutadiene. 
Whereas Smith and co-workers in 1999 were still limited to 100 carbon atom simulations 
in 2005 Tsolou and co-workers were able to study the self-diffusion and dynamic 



32 



structure factor of polybutadienes with up to 400 carbon atoms |33). 

Despite the fast growth in computing power, the extremely fast increase of the 
longest relaxation time with chain length, n^'^, implies that atomistically detailed MD 
simulations will remain limited to only slightly entangled samples. Longer chains can 
and have been studied, for example Ryckaert studied a melt of fully atomistic Ciooo 



chains 34 , but this can only be done for simulation times much shorter than the longest 



relaxation time. 



3. Coarse-graining: lumping together the atoms 

3.1. Theory of coarse- graining 

In order to increase the time and length scales accessible in the simulation of polymers, 
detailed atomistic models are replaced by coarse-grained models in which each particle 
represents a collection of atomic particles. In this review we focus on the bottom-up 
{ah initio) approach where the interactions between the coarse-grained particles are 
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directly linked to the atomistic interactions and are aimed at correctly reproducing the 
structural, thermodynamic and/or dynamical properties [3|[7|[9| [l0|[35] . 

When coarse-graining polymers, the polymer chain is subdivided into subchains, 
partitioning the degrees of freedom in two sets: the coarse-grain coordinates and 
momenta {i?, P}, which are the centres of mass (or characteristic atom) positions and 
momenta of the subchains, and the 'bath' coordinates and momenta {g, p}, which are the 
remaining internal coordinates and momenta describing the details of the configurations. 
The potential of mean force governing the forces among the coarse-grained coordinates 
is given by 

A{R) = -ksT In j dp dg exp {-(3Hb{p, g; R)) , (1) 

where /3 = 1/ (ksT) and Hb{p, g; R) is the bath Hamiltonian, equal to the sum of kinetic 
energy Tsip^q) of the bath variables and the potential energy ^{R,q) of the entire 
system. In a system containing only coarse-grain particles, employing the potential 
of mean force A{R) ensures a correct distribution of R coordinates, as well as correct 
thermodynamic properties. Unfortunately, the potential of mean force is generally a 
complicated function of all coordinates R, and invariably includes complicated multi- 
body interactions. For practical and computational reasons it is impossible to calculate 
and store the potential of mean force for all possible multibody configurations. Some 
approximations need to be made. The most widely used approximation is that the 
potential of mean force can be represented by pairwise and sometimes triplet terms. 
We would like to issue a warning that a naive use of the interactions at the coarse- 
grained level can lead to incorrect values for some thermodynamic properties. More 
specifically, the thermodynamics of the coarse-grained system will not correspond to 
the thermodynamics of the microscopic system if one falsely assumes that one is still 
dealing with an atomistic system, where the coarse-grain particles are the only particles 
present and the effective pair interactions are the only energy terms present. It should 
be intuitively clear that the relatively softer interactions will then lead to a pressure 
which is generally too small and a system which is too compressible. The source of this 
thermodynamic inconsistency can be traced to the state point dependency as well as the 
pair approximation of the potential of mean force. For a discussion on these matters. 



the reader is referred to 36 38 



The theory underlying coarse-graining of dynamics is the Mori-Zwanzig formalism 



10 ,35 ,39 ,40 , where a projection operator technique leads to a Generalized Langevin 
Equation governing the dynamics of the coarse variables. Although the dynamics of 
the atomistic model is deterministic and conservative, the dynamics of the Generalized 
Langevin Equation is stochastic and includes dissipation and memory effects: 

-^it) = -M-^t) + -T^ldr PUt - r)Un (r; R{t - r)) . (2) 

Here P„ is the momentum of coarse-grained particle n, A is the above-mentioned 
potential of mean force, is a random force, and Cmn is a time-dependent friction 
matrix. The Mori-Zwanzig formalism provides an exact derivation for this friction 
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matrix, but its interpretation is not straightforward: Cmn (t; R{t — r)) is the correlation, 
at the point where the coordinates R equal R{t — r), of the random force at time zero 
on particle m with the random force propagated by a complicated operator on particle 
n. Specifically, this operator is exp{(l — V)iCT}, where C is the Liouville operator and 
V projects onto the coarse-grain coordinates [Toj. If the theory is followed to the rule 
(including the exact potential of mean force), a correct description of the structure and 
dynamics follows automatically. However, as the above equation shows, the resulting 
friction on a particle again depends on all coordinates R and momenta P. Compared to 
the case of the potential of mean force the situation is aggravated, because the friction 
depends also on the coordinates and momenta in the past. Clearly, for practical and 
computational reasons approximations need to be made. For example it is often assumed 
that the friction is pairwise additive with the contribution of a pair depending only on 
the positions and momenta of the two particles involved. Also memory effects are often 
ignored. This may or may not be correct. Making the right approximations, often based 
on physical intuition, is the responsibility of the coarse-grainer. 

3.2. Classification of coarse-grained dynamics methods 

The term 'coarse-grained dynamics' is used rather loosely in the literature. Many 
different levels of coarse-graining are encountered, with different approaches to treating 
the friction. In our view it is important to distinguish between two main classes, 
which may be termed coarse-grained molecular dynamics and coarse-grained stochastic 
dynamics, see Fig. [T| The difference is the following: 

• Coarse-grained molecular dynamics (CGMD) applies to simulations where 
a few atoms (up to one chemical repeat unit or five backbone carbon atoms) 
are coarse-grained to relatively hard beads. In this case the additional friction 
components are often not dominant, and in most cases are ignored altogether. 
Ignoring the friction leads to too fast dynamics, which is beneficial for equilibration, 
but requires a rescaling of time when quantitative comparison with real dynamics 
and rheology is desired. 

• Coarse-grained stochastic dynamics (CGSD) applies to simulations where 
many atoms (many chemical repeat units or more than 10 backbone carbon atoms) 
are coarse-grained to relatively soft beads or 'blobs'. The friction and stochastic 
forces dominate the interactions, and cannot simply be ignored but must be included 
in the equations of motion. Often it is assumed that the stochastic forces are 
delta-correlated, i.e. have no memory. Unless countermeasures are taken, the 
soft interactions cannot prevent crossing of the bonds between the coarse-grained 
particles, leading to unrealistic dynamics for entangled polymers. 

In the following sections we will treat both levels of coarse-graining. 




Figure 3. (color online) Coarse-graining of polystyrene to two types of beads 'A' and 
'B' representing part of the backbone and a phenyl-ring, respectively. Dashed lines 
show the bonds between the coarse-grained beads A and B. Reproduced from j4f]. 



4. Lumping a few atoms into beads 

Let us first focus on coarse-grained molecular dynamics (CGMD) simulations. A coarse- 
grained particle may represent the centre-of-mass of a small group of atoms, or the 
position of a characteristic atom which is important in determining the local structure 
of the polymer backbone or a sidegroup. An example is given in Fig. |3]where polystyrene 
chains are coarse-grained to 'A' and 'B' beads representing a part of the backbone and 



the phenyl-ring, respectively 41 



4.I. Targeting the structure of real polymer chains 

Most coarse-graining efforts focus at reproducing the structure of the underlying 



microscopic chain 3 5[]8 42 



The structure of a polymer melt can be described in terms of distribution functions 
of distances and angles between the coarse-grained particles. Once it is decided what 
the position of a coarse-grained particle represents, it is rather straightforward to 
measure distribution functions for these coarse-grained coordinates in atomistically 
detailed molecular dynamics simulations. Intramolecular distribution functions include 
those for the distance between two adjacent coarse-grained beads, the angle between 
three consecutive beads, and the dihedral angle between four consecutive beads. Other 
structural distribution functions include those for distances between beads on different 
chains, and can also include those for distances between beads on the same chain but 
sufficiently far away along the backbone. These distribution functions are used as target 
functions for the coarse-grained interactions. 

Usually, the local intramolecular interactions between consecutive beads are 
relatively strong and therefore rather unperturbed by intermolecular interactions. In 
such simple Boltzmann inversion is enough to obtain the coarse-grained 

interaction [5]. For example, from a distribution Pbond(^) of bond lengths r between 
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Figure 4. Nonbonded A-A pair distribution function for a polystyrene melt (see 
Fig. [3]) at T = 463 K, obtained from atomistic molecular dynamics simulations (solid 



line) and coarse-grained simulation (dashed line). Reproduced from 44 



consecutive coarse-grained particles (obtained in an atomistic simulation and already 
corrected for the Jacobian 47rr^), the coarse-grained bond potential is derived as 

Vhond{r) = -ksT In Pi,ond{r). (3) 

One word of caution is at place here. Factorising the various distribution functions into 
independent parts representing bond distances, bond angles and dihedral angles ignores 



any correlations which may be present between them 43 , 44 



The nonbonded interaction, represented by the radial distribution function g{r) 
between nonbonded coarse-grained beads, is however much more subtle. It contains both 
entropic and enthalpic contributions and a simple Boltzmann inversion is insufficient 
because this ignores important packing effects [4|[5||42||45]. Different methods have 
been proposed to find the nonbonded interaction which reproduces the target radial 
distribution function gta,rget{r)- One option is to choose a simple functional forms for 
the interaction, based on physical intuition, and optimise the parameters of this function 



to reproduce the complicated radial distribution function as close as possible 718,44,46 



For example, the nonbonded interactions between the coarse-grained polystyrene beads 
of Fig. |3] have been found by optimising the parameters of a generalised Lennard- Jones- 



type interaction 44 . As Fig. |4| shows, the resulting nonbonded radial distribution 
functions agree very well with those obtained from atomistic simulations. 

In some cases a predefined functional form is not fiexible enough to capture the 
intricacies of the target radial distribution function. More complicated analytical 
potentials, defined piecewise for different ranges, allow more accurate reproduction of 
the complicated distribution functions 46 48 . The parameters of these models may be 
obtained semi-automatic or fully-automatic using simplex optimisation [5 ,49 52 . 

Another popular technique to obtain the coarse-grained interactions is Iterative 
Boltzmann Inversion (IBI) [5||47||53|[54]. Focusing on the non-bonded interactions, for 
every iteration a one-to-one correspondence between the effects at a distance tq and the 
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non-bonded potential V{rQ) at the same distance is assumed. In the beginning, a starting 
potential Vo(r) has to be guessed, which is usually the Boltzmann inversion of the target 
distribution function. At each iteration stage n, a simulation is performed with the last 
guess Vnir) of the non-bonded potential, and the radial distribution function gnij) is 
measured. This leads to a new guess for the non-bonded potential: 

Vn+i{r) = V;(r) + ksTln (4) 

5'targetV j 

The iterations continue until the difference between gn{'r) and (^target is deemed 
satisfactory. Usually only a few iterations are necessary [5|[47]. 

Other techniques have been used as well. Within our own group we have 
decribed coarse-grained interactions by a number of functions with coefficients which 
are left for automatic optimisation by treating them as dynamical parameters in an 



extended ensemble 55 . Notably, Ashbaugh and co-workers used a combination of 
molecular dynamics simulations and inverse Monte Carlo methods to self-consistently 
map structural correlations from atomistic simulations of alkane oligomers onto coarse- 



grained potentials 56 



A warning is at place here. All techniques described above aim at reproducing 
the structure by means of pair interactions. This means that the true potential of 
mean force, Eq. ([T]), is only approximated and many-body interactions beyond a certain 
number of particles (usually two) are ignored. As a consequence, it is not guaranteed 
that the thermodynamic state is correctly described [5 ,47,55 . 



An option is to include thermodynamic properties in the optimisation scheme. For 
example Reith et al. 47 added an extra pair-potential with an absolute value decreasing 
linearly to zero at a cut-off range. The amplitude (positive or negative) and cut-off range 
were used as new parameters for re-optimisation. This yielded a radial distribution 
function which did not deteriorate too strongly, combined with a correct pressure. Also, 
Nielsen et al. 57 have developed a coarse-grained model for n-alkanes (polyethylene) by 
optimizing bond and bend parameters and the non-bonded Lennard- Jones parameters 
to match observables from atomistic simulations as well as experimental surface tension 
and bulk density data. 

As already mentioned, when coarse-graining only a few atoms to beads the 
interactions are still relatively strong. This means that, even though longer simulation 
times can be reached, equilibration is still sometimes an issue. In such cases the 
equilibration techniques developed for atomistically detailed models can be used with 
high efficiency. For example, recently End-Bridging techniques have been applied to 



coarse models of poly (ethylene terephtalate) and atactic polystyrene 23-25,58 



In summary, when coarse-graining one should always keep an eye on what properties 
one is interested in. Because in practice one is limited to approximating the potential 
of mean force in terms of pair interactions, it is impossible to correctly represent both 
the structural and thermodynamic properties. Improving the agreement with one will 
invariably deteriorate the agreement with the other, so compromises need to be made. 
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Jf..2. Lattice models 

Up to this point we have only considered molecular models in continuous space. 
Atomistic models have also been mapped to lattice models. One of the most important 
lattice models for polymers is the Monte Carlo-based Bond Fluctuation (BF) model 



59 61 , where each super-atom occupied eight sites on a cubic lattice. Bond distances 
and angles can vary between different discrete values. The advantage of this approach 
is that the possible distances and angles can be chosen in such a way that bond crossing 
becomes impossible. Mapping of real polymers onto the BF model has been tried for 



bisphenol polycarbonates and polyethylene 62 65 . Coarse-graining using the BF model 
is extensively reviewed in [s]. 

Mattice and coworkers have developed a mapping from a rotational isomeric 
state (RIS) description of a polymer onto a diamond lattice for the second nearest 



neighbour positions of polyolefine backbones 66 69 . The intermolecular potential in 
this model, with different values of the energy for different ranges of distances between 
the superatoms, was tuned by aiming at producing a zero second virial coefficient 
(corresponding to theta-conditions) as well as a correct radius of gyration for polymers 
such as polyethylene and polypropylene. For the dynamics local Monte Carlo moves, 
including a crankshaft move, were introduced. Recently the dynamical updates have 



been redesigned to study the mean-square-displacement in polyethylene melts 70 and 



in bidisperse polyethylene melts 71 



4.3. Backmapping: reinserting atomistic details 



An important task of CGMD models is to quickly generate well-equilibrated atomistic 
polymer structures. For this it is of course necessary to be able to map an equilibrated 
coarse-grained model back onto an atomistic model. Such a backmapping procedure is 
possible, as has been extensively reviewed in [3|[7|. 

Among the first polymers to be equilibrated in this way are various polycarbonates 
46 , 48 , 72 , 73 , but more recently backmapping procedures have also been applied to 



polystyrene 44 ,58 ,74 -76 and polyamide 77,78 chains. 



An exciting new development, with possible applications to determination of 

In 



the rheology, is the extension of backmapping to non-equilibrium situations 79 



this method deformed conformations are maintained during backmapping by applying 
position restraints. The method has been demonstrated for atactic polystyrene under 
steady shear fiow. 



4.4- Dynamic properties from coarse-grained molecular dynamics: rescaling of time 

CGMD simulations have not only been used to accelerate the equilibration of polymer 
melts, but also to study their dynamics (7|[9]. We have already discussed that additional 
friction factors and stochastic forces are inexorably linked to the procedure of coarse- 
graining. However when the degree of coarse-graining is not too large, friction is often 



Coarse- graining of the dynamics of entangled polymer melts 



13 



ignored, leading to too fast dynamics. This may optimistically be referred to as 'speed- 
up'. In such cases, when dynamic properties are presented, it is assumed that all 
important processes are accelerated by the same factor (usually 2 to 5), i.e. that a simple 
rescaling of time by a single scaling factor can recover the real dynamics. Although in 
most cases this assumption has no theoretical justification, it nevertheless seems to be 
correct for low degrees of coarse-graining, as may be witnessed by the many successful 



applications of time rescaling that may be found in the literature 41 , 44 53 , 74 , 75 , 80 - 88 

Calibrating the time-scale may be done by comparing chain diffusion coefficients, 
although for most realistic polymers the atomistic simulations do not reach long enough 
times to observe free chain diffusion. Since the coarse-graining level is relatively modest, 
there is hardly any difference between a bead in a relatively short (but still polymeric) 
chain and a bead in a very long chain. Thus it can be expected that a diffusion-based 
calibration may be done using relatively short chains. Harmandaris and co-workers have 
been very successful in consistently determining the time scaling factors by comparing 
the atomistic and coarse-grained mean-square-displacements of relatively short chains 



of polystyrene |41||44|[74] . Fig. [5j shows a recent example 41 , where it is found that the 
mean-square displacement of the chain centre-of-mass and the end-to-end vector time 
correlation can both be accurately reproduced by the same scaling factor S. It is also 
possible to use other shorter local timescales to calibrate the time-scale. For example 
one can use the time associated with the relaxation of subchains, as sampled by time 
correlations of the higher Rouse modes [l]. 

5. Lumping many atoms into blobs: the role of friction and uncrossability 

In order to reach even larger time and length scales one may combine many more atoms, 
say > 10 monomeric units, into one coarse-grained particle. The particle positions must 
be updated according to a coarse-grained stochastic dynamics (CGSD) scheme because 
the friction and stochastic forces often dominate the interactions. Neglect of friction and 



stochastic forces can lead to erroneous dynamics 89 ,90 . For example, the motion at 
certain length scales may appear as oscillatory, whereas the real polymer is overdamped 
at these scales. 

5.1. Friction and stochastic forces 

As already mentioned, the Mori-Zwanzig formalism of coarse-graining the dynamics 
leads to a Generalized Langevin Equation (CLE), Eq. (|2|, in which the friction and 



stochastic forces appear with memory 10 39 40 91 93 . Even though this equation 



is exact, it is very difficult to use in practice. In the first place, evaluation of the 
memory terms requires taking averages over a projected dynamics which one does not 
know exactly how to generate. In the second place, integrating the resulting integro- 
differential equations is an extremely challenging numerical task. 

Usually one therefore assumes that the time scale associated with changes in the 
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time (ns) 

Figure 5. Example of a successful time mapping of coarse-grained (CG) simulations 
onto chemically realistic united-atom (UA) molecular dynamics simulations of 
polystyrene. Two melts, with molecular weights 1 and 2 kDa have been studied, 
both at T = 463 K. The time scaling factor S needed for the chain centre-of-mass 
mean-square-displacement (a) is consistent with the time scaling factor for the time 
correlation of the end-to-end vector (b). Reproduced from |41| . 



positions of the coarse-grained coordinates is much larger than the time scale associated 
with the decay of the friction memory. In this so-called Markov approximation the time 
dependence of the friction in Eq. ^ is represented by a delta-function, resulting in 
an ordinary Langevin equation for the dynamics of the coarse-grained particles. Whether 
or not the Markov approximation is valid must be carefully evaluated in each particular 
case. 

Another difficulty remains: the friction and stochastic forces on a particle still 
depend on the positions and velocities of all particles. As already mentioned, it is often 
assumed that the friction is pairwise additive with the contribution of a pair depending 
only on the distance and relative velocities between the two particles involved. Doing 
so, we arrive at the pair friction model employed in the so-called dissipative particle 
dynamics (DPD) method 94 . Although this method is still relatively costly, especially 
when the frictions are high, it is rather popular since it conserves momentum, which 
is a prerequisite for hydrodynamic behaviour, and it is easy to use in non-equilibrium 
simulations. Within our group we have studied in detail the pair friction parallel and 
perpendicular to the axis connecting the centres-of-mass of two halves of a polymer chain 



in a melt 35 . Also, coarse-grained forcefields for polyethylene and cis-polybutadiene 
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for use in DPD have been determined in the group of Rousseau 95 96 . The Rousseau 



group discussed the dependence of the friction coefficient on the coarse-graining level in 
view of the overall scaling of the dynamical properties. 

The use of pair frictions in DPD is necessarily linked to the use of a second order 
propagator. This implies that one must use low frictions to reach long time scales 
(high frictions require too small time steps). The damage an artificially lowered friction 
does to the rheological properties is often counterbalanced by using harder conservative 
interactions [97^|. This often leads to unrealistic structure on the scale of the bead. 

Fortunately, in the case of polymer melts it is often admissible to use a scalar 
friction (i.e. with a static background), because the friction may be thought of as being 
caused by the motion of a (part of a) chain relative to the rest of the material, which 
to a first approximation may be taken to be at rest. Propagation of a velocity field 
as in a normal liquid is highly improbable, meaning that hydrodynamic interactions 
are screened [l]. The positions of the coarse-grain particles can then be updated by 
a first order Brownian dynamics propagator, using much larger time steps than are 
possible with a second order propagator. One should be aware, however, that such a 
first order Brownian dynamics propagator is not Galilean invariant. Extra care should 
therefore be taken when simulating flow of a polymer melt. Often it is assumed that 
the particles flow affinely with the imposed deformation, although in the case of shear 



flow this assumption may be relaxed 98 -100 



5.2. Uncrossability of bonds and entanglement 

The advantage of lumping together only a few atoms in CGMD models is that the 
resulting interactions are still sufficiently strong to prevent bond crossing. In bead- 



spring models, such as the important FENE model 101 111 , uncrossability is ensured 



by the use of relatively hard beads. Such models therefore belong to the class of CGMD 
models. Relatively small friction and stochastic forces have been added to bead-spring 
models, but with the goal of stabilising the integration of the equations of motion over 



very long time scales 102 , rather than representing the physical friction caused by the 



degrees of freedom that have been coarse-grained out. 

As more and more atoms are coarse-grained into one particle, the interactions 
between them become softer and softer 88 , 112 , 113 , see Fig. [61 Beyond a certain degree 



of coarse-graining, bonds will be able to cross each other. Without additional measures 
the important entanglement effect, leading to altered and much slower dynamics, is 
therefore lost. 

An algorithm that can detect and prevent unphysical bond crossings has been 
developed by us 112 . Employing this so-called Twentanglement algorithm to a 



simulation of highly coarse-grained polymer chains will reintroduce the entanglement 
effect. The principle of the algorithm is depicted in Fig. [7} The bonds are considered 
to be elastic bands between the bonded particles. The algorithm keeps track of all 
(unattached) bond vectors which are close together. For each bond vector and at each 
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Figure 6. The potential of mean force between bonded (circles) and nonbonded 
(squares) coarse-grained pieces of polyethylene (T = 450 K), each piece representing 
the centre-of-mass of 20 carbon groups. These potentials were obtained from 
distribution functions measured in atomistically detailed molecular dynamics 
simulations of Ci2oH242- Note that fc^T = 3.74 kJ/mol, meaning that without 
additional measures bonds can easily cross each other. 




Figure 7. (color online) Principle of the Twentanglement algorithm, (a) Two line 
segments representing a bond are closing in on each other. (6) At a certain moment 
these bonds will touch. An 'entanglement' is created at the crossing point X. (c) After 
this, the bonds are viewed as slippery elastic bands. The elasticity will slow down the 
relative speed of the bonds. This sequence of events may also be reversed. 



instant of time the following triple product is calculated: 

Vij = (ri - Vj) ■ [(ri+i - Ti) X (r^+i - rj)] , (5) 

where we have used Fig. [T] as a reference for the indices. The absolute value of Eq. (jsj) is 
the volume of the parallelepiped defined by the vectors rj+i j, r^+i and Aside from 
some pathological cases 112 , if Vij changes sign from one time step to the next, a bond 
crossing may have occurred. Additional checks are made to ensure that the crossing is 
taking place along the physical part of the line segments (the above equation checks if 
two infinite lines have crossed). If a real bond crossing has occurred, an entanglement 
is created at the crossing point. Subsequently, the associated volume Vij will serve to 
detect future disentanglements. If the volume Vij of the four objects surrounding an 
entanglement changes sign, a possible disentanglement has occurred, i.e. Fig. [7] may 
also be read backwards. Usually only a few of the uncrossability constraints contribute 
to an entanglement in the usual sense of longlasting obstacles, slowing down the chain 
movement. For instance, a C60H122 chain is generally considered not to be entangled. 
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Figure 8. The zero-shear viscosity ryo as a function of molecular weight for 
polyethylene at T = 450 K from coarse-grained simulations using Twentanglement 



114| (solid squares) and experiment [119| (dashed lines) . The result from a chemically 

is also included (open circle). Observe that the steep 
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detailed simulation of C 
increase in viscosity does not start at the entanglement molecular weight Mg but at a 



critical molecular weight Mc which is several times higher 114 



yet many 'entanglements' occur in a coarse-grained simulation. The dynamics of the 



'entanglement' points is determined by a balance of forces, as described in detail in 112 



Because uncrossability is explicitly taken care of, there is much freedom to choose 
the number of monomers per coarse-grained particle. In our previous work, the following 
considerations were taken into account: (1) the degree of coarse-graining should be 
large enough to allow for a significant increase in the time and length scales accessible 
to simulation, and (2) if one wants to study reptation effects, the degree of coarse- 
graining should not be so large that the size of a coarse-grained particle exceeds the 
typical diameter of the tube in the reptation picture, in other words the entanglement 
length. A suitable choice for polyethylene was to coarse-grain 20 CH2 units to one 
particle 112 . We have also been able to measure the required friction on such a coarse- 



grained particle from a chemically realistic simulation by analysing the correlations 
in the constraint force required to keep the location of a coarse-grained particle fixed 
in an atomistically detailed simulation 114 . By using the Twentanglement algoritm 



we have been able to predict the dynamics and rheological properties of polyethylene 
melts up to Ciooo (Fig. [s]) 114 117 and the dynamics of poly(ethylene-alt-propylene) 
melts 118 . We also applied the method to entangled wormlike micelles of a specific 



surfactant chemistry, and found quantitative agreement with experimentally determined 



rheology 98 



Analysing the results of these coarse-grained simulations we found that the chain 
stiffness is an important ingredient for the dynamics of polymer chains. Indeed, Faller 
and Miiller-Plathe noted that chain stiffness intensifies the reptation characteristics of 



polymer melt dynamics 81 , 120 
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Topological constraints have also recently been introduced in DPD simulation 
In standard DPD, the conservative forces between the polymer 



methods 97,113 121 



segments decrease linearly with increasing pair distance. The interactions are therefore 
naturally quite soft and polymer chains built from such soft beads behave as phantom 
chains who pass freely through each other 122 . Nikunen et al. prevented such chain 



crossing by simply increasing the strength of the DPD forces 97 , whereas Liu et al. 



added a rigid core to the DPD spheres 121 . Note that both these schemes essentially 



make the conservative interactions more rigid, which is opposite to the trend observed 
when coarse-graining to higher levels. One of the dire consequences is the introduction 
of unwanted structure on the scale of the DPD particle, which makes these methods 
unsuitable for high degrees of coarse-graining. 

Uncrossability of chain bonds can alternatively be introduced by introduction of 
an additional repulsive interaction which is based on the distance of closest approach 
between two bonds 113 ,123 ,124 . Another possibility to conserve the topology, at least 
for Brownian dynamics simulations, is to simply forbid random displacements that lead 
to crossing of bonds 125 . In a similar vein, one of us has very recently developed 



an algorithm which efficienly treats collisions between infinitely thin bonds in such a 
way that momentum and energy are conserved locally 126 . This allows the study 



of the relative importance of topological and hydrodynamic interactions. It should 
be noted that in all methods described in this paragraph the bonds themselves stay 
rigid. A collision between two such bonds is a hard collision. In real polymer melts the 
bonds between consecutive coarse-grain beads are highly flexible, and the collisions will 
consequently be much softer (e.g. as in Twentanglement). As noted above 81,120 , the 
intrinsic stiffness of the bonds may lead to artifically enhanced entanglement effects. 
These methods are therefore more suitable for simulation of semifiexible chains, or for 
flexible chains at such a low degree of coarse-graining that local chain stiffness is still 
important. 



6. Primitive paths 

We have seen in the previous section that beyond a certain level of coarse-graining it 
becomes necessary to somehow maintain the information about the topology of the 
system, otherwise the important entanglement effect is lost. One option is to use 
algorithms such as Twentanglement. Another option is to switch to a one-chain model 
based on the reptation model [l]. This however poses several problems. Reptation 
theory is still facing many questions which cannot be answered within the model. 
For example, to quote a recent review of Likhtman p], it is still not clear what is 
the entanglement or tube field restricting the chain motion. What are the statistical 
properties of the tube - is it semifiexible on the length scale of the tube diameter? What 
is the microscopic basis for the assumptions used to describe constraint release in the 
linear and non-linear regimes? What happens to the entanglements and the tubes after 
large deformations or in fast flows? A practical question of relevance for quantitative 
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predictions is how the numerical parameters of the tube model, such as the entanglement 
mass and tube diameter, can be related to the chemical composition of the polymer. 

6.1. Primitive path analysis 

The questions posed above have provoked an outburst of activities which try to 
determine the primitive path and entanglement characteristics from more detailed 
simulations [6|[2^[M|[58|[8il [l02|p7Hl34l ■ AH efforts are based on Edward's view [I] of 



a primitive path as the shortest path remaining when one holds chain ends fixed, while 
continuously reducing (shrinking) a chain's contour without violating uncrossability 



130 . The chain contours are reduced simultaneously until topological constraints block 
further shrinkage. 

Everaers et al. introduced one such minimization method to extract the tube 
diameter and entanglement length from more detailed simulations of linear chains 



127, 128 , see Fig. |9j On the basis of the tube model, both the real chain and the 
primitive path beyond some length scale are described as random walk chains with 
the same mean square end-to-end distance (-R^), but with different contour lengths. 
The entanglement length and tube diameter d correspond to the length and end-to-end 
distance of an entanglement strand, which is identified with the Kuhn segment of the 
primitive path. Thus it is assumed that the following relations hold (!]: 

<i=^. 'V.-ivi^. (6) 

where (Lpp) is the average primitive path contour length and is the number of Kuhn 
segments of the real chain. Nf, is the entanglement molecular weight expressed as the 
number of Kuhn segments of the real chain. A refinement of this method to exclude 



self-entanglement was introduced later 135 . It turned out that entanglements between 



distant sections of the same chain make a negligible contributions to the tube and that 
the contour length between local self-knots is significantly larger than the entanglement 
length. 



The above method has been applied, among others, to polystyrene melts 41 



to semidilute solutions of stiffer chains such as actin 136 , and to bisphenol-A 



polycarbonate 82 , which was found to have a surprisingly low entanglement length. 
The estimated entanglement molecular weight Me is usually in good agreement with 
the experimental value obtained from the plateau modulus Gq. 

The plateau modulus is a mechanical property of the system. Following a small step 
strain 70, the shear stress 7oG(t) initially decays according to local relaxation processes 
on length scales smaller than the tube diameter. For well-entangled chains, after the 
entanglement time Te, further relaxation of G{t) is severely delayed by the topological 
interactions, leading to an apparent plateau; only for times much longer than can 
G{t) finally relax to zero. The plateau modulus may be defined as the modulus at 
times larger than Tg, but much smaller than the terminal relaxation time. The relation 
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Figure 9. (color online) Result of a primitive path analysis of a melt of 200 chains each 
consisting of 350 beads. The primitive path of one chain is shown (in red), together 
with all of those it is entangled with (in blue). The primitive paths of all other chains 
are shown as thin lines. Reproduced from [127|. 



between the plateau modulus and the entanglement molecular weight is given by [T]: 

Here p is the polymer (mass) density. This expression is very similar (apart from the 
factor 4/5) to a prediction for the affine model of rubber elasticity [l]. Indeed, the good 
agreement between the entanglement length from the primitive path analysis and the 
entanglement length from the plateau modulus may intuitively be understood because 
they both probe the same (albeit temporary) rubber-like elasticity of the underlying 
primitive network of chains. 

The agreement with other estimates of the entanglement length is less good 
however, especially when these estimates are based on characteristic changes in the 
time- dependence of quantities which measure the motion of individual chains within 
their prospective tubes 61,116j. It may be argued that this disagreement is caused 



by the dependence on non-exact theory in interpreting these measurements 111,130 



This is certainly a factor, but the remarkable observation is that with reasonable 
assumptions in most cases time-resolved measurements yield an entanglement length 
which is consistently larger than the one based on the plateau modulus. For example 
for polyethylene, entanglement length estimates from the dynamic structure factor, the 
cross-over in Rouse mode time correlations, and the entanglement time Te (estimated 
directly from the time of the inflection point in G{t)) all give a value which is consistently 



larger by a factor of 1.5 compared to the estimate from the plateau modulus 116 



The disagreement with the 'plateau modulus' entanglement length may have two 
causes. First, dynamic time- resolved quantities are probably more sensitive to details of 
the statistics (e.g. fluctuations) of the entanglements. Second, dynamic time-resolved 
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quantities may be more sensitive to the finite stiffness of a real polymer chain on the 
scale of the tube diameter. 
In a recent article 



137 



Greco noted that, because the tube diameter is 
relatively small, it is very important to include fluctuation effects in the description 
of entanglements. Unfortunately, at the time of writing this review, there is still 
no definite agreement on the best way to determine the primitive paths and their 
statistics. Different coarse-graining schemes for arriving at the primitive path lead to 
differences in the statistics. For example, Kroger and co-workers proposed an alternative 
algorithm which returns a shortest path and related number of entanglements for a 



given configuration of polymers 129,133 . Primitive paths are treated as infinitely thin 



and tensionless lines rather than multibead chains, and excluded volume is taken into 
account without a force law. Their implementation allows construction of an optimal 
shortest path for 3D systems. The number of entanglements is then obtained from 
the shortest path as either the number of interior kinks, or from the average length of 
a line segment. This primitive path analysis has been applied to polyethylene C24 to 
Ciooo 133,138,139 , as well as the FENE model 133 , vahdating analytical predictions 



of Schieber 140 about the shape of the distribution for the number of strands in a 



chain at equilibrium. It was also concluded that the number of entanglements obtained 
by assuming random walk statistics, as is done in the work of Everaers and co-workers, 
deviates significantly from the predictions by the algorithm of Kroger and co-workers. 

The random walk approximation partly fails because real polymers are typically 
still somewhat stiff at the scale of the average distance between two uncrossability 

As a result, the tube diameter cannot be equated to the Kuhn 



130 133 



constramts 

step length of the primitive path. Although this fact is well-known for semifiexible 
fiber suspensions, where theoretical expressions for the tube diameter are different from 
expressions for the deflection length (average distance between successive collisions of the 
fiber with its tube) [141 , it seems to be often ignored in the polymer melt literature. 
Stiffness effects have been introduced recently in primitive path analyses by Hoy et 



al. 133 and by Tzoumanekas et al. 132 



6.2. Primitive chain networks and slip-link models 

Although the precise characteristics of entanglements and/or confining tubes have 
remained somewhat elusive, this has not stopped researchers from trying to predict 
polymer melt rheology by employing primitive path-based models. 

An important development in this direction is the primitive chain (or slip-link) 
network model developed by Masubuchi and others 142 , 143 . In this model polymer 



chains are coarse-grained to the level of segments between entanglements, see Fig. [TO] 
Chain coupling is achieved by confining two chains by a common slip-link placed at a 
specific position in 3D space. The chains fulfill a force balance at the slip-links, similar 



to the treatment of uncrossability constraints in the Twentanglement algorithm 112 



which act on a somewhat smaller scale. The Langevin equation for the slip-links contains 
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Figure 10. (color online) Schematic representation of a primitive chain network 
model, (a) The motion of slip-links (rings) is influenced by the tension in the chain 
segments (arrows) and an osmotic force, (b) The motion of the monomers through 
the slip-links (arrows) results in reptative motion. (Picture kindly provided by Y. 
Masubuchi.) 



both the tension in the chain segments emanating from the shp-hnk and an osmotic 
force arising from density fluctuations. The motion of monomers through the shp-hnks 
ultimately generates reptation as well as tube length fluctuations. When creating new 
slip-links (i.e. new entanglements), a partner chain segment is chosen randomly among 
those spatially close to the advancing chain end. 

The primitive chain network model correctly predicts the scaling of the longest 
relaxation time and the self-diffusion with chain length 142 . The model was later 



modified to also include the concept of dynamic tube dilution 144 , and good agreement 
between simulated and real linear and nonlinear rheology was found. Extensions to 
branched polymers, blends and block copolymers have all been made 145 146 . Recently 



the model has been subjected to various sudden deformations to study the damping 
function, which is essentially the ratio of nonlinear stress relaxation to the linear shear 
relaxation modulus 



147 149 



It was found that the force balance is a dominant 
correction over the basic Doi-Edwards theory as compared with the effect of convective 
constraint release. Furthermore, the predicted normal stress ratio in shear, a quantity 
which is very sensitive to different modelling assumptions, was found to be in good 
agreement with experiments. 

In order to facilitate fast computations of the rheology, various single-chain slip-link 
models have also been developed 150 155 . In such models an ensemble of primitive 
paths is followed in time. Coupling between chains is not explicit, as in the above 
primitive chain network model, but virtual by defining random relations between slip- 
links on different chains. Whenever a slip-link of a chain is destroyed due to tube 
renewal at a chain end, the corresponding partner along the contour of another chain is 
also eliminated. The same applies to creation of new slip-links. 

Very recently, the effectiveness of such a single-chain slip-link model in describing 
the dynamics and rheology of entangled polymers was tested by systematically 
comparing slip-link with FENE molecular dynamics results 156 . The parameters of the 
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model were determined by using the mean-square displacements of one particular chain 
length as a target function. Using the same set of parameters, it was then tested if the 
predictions of the mean-square displacements for other chain lengths agreed with the 
MD calculations. This was followed by a comparison of the shear relaxation moduli G{t). 
The work identified a limitation of the original slip-spring model in describing the static 
structure of the polymer chain as seen in MD, which was remedied by introducing a 
pairwise repulsive potential between the monomers in the chains. Also, poor agreement 
of the mean-square monomer displacements at short times could be rectified by the use 
of generalized Langevin equations for the dynamics. The latter of course somewhat 
diminishes the computational efficiency of the method. 

Slip-link models are generally successful in their overall 'rheological' performance 



and predictions 130 . They also offer the possibility to introduce new tube renewal 
mechanisms and, by turning these on or off at will, probe the effect on the rheology 
independent of other mechanisms. A disadvantage of both single-chain and multi-chain 
slip-link models is that the link with the chemistry usually is far away; if a comparison 
with a microscopically detailed model is made at all, the parameters are usually obtained 
by fitting rather than an ab-initio coarse-graining procedure. 

An attempt at systematically bringing in structure from a detailed model, albeit 



again a bead-spring type model, was introduced by Rakshit and Picu 157 . In their 
method reptation is enforced because the chain inner blobs are constrained to move 
along the backbone of the coarse grained chain (the primitive path), while the end blobs 
move in the three-dimensional embedding space. These end blobs continuously redefine 
the diffusion path for the inner blobs. It was shown that various distribution functions, 
relaxation rates, and the diffusion dynamics are properly represented. Recently, this 



model was used to study start-up and step strain shear flows 158 . The authors 
showed that their method reproduces several features observed experimentally such as 
the overshoot during start-up shear flow and shear thinning at large shear rates. 

7. Super coarse-graining: a polymer as a single soft particle 

Industrial processing of polymers usually involves extrusion of a hot melt through a 
die or injection of a melt into a mould. Modern applications also include nanoparticles 
in the melt to form composites. To optimise the process parameters for such complex 
flows it is important to understand the interplay between time and length scales of 
the polymer and those imposed by the complex geometries. The low frequency linear 
and non-linear rheology of polymer melts is basically concerned with displacements of 
centres of mass only. Given the need to also reach large length scales, it is therefore just 
natural to try to develop simulation models that only keep track of the dynamics of the 
centres of mass of all diffusing constituents. 



An early application of this idea was presented by Murat and Kremer 159 . These 
authors developed an efficient and rather general model in which whole polymer chains 
are represented as soft ellipsoidal particles. The interaction between two such particles 
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is taken to be proportional to the spatial overlap of their monomer density distributions. 
Since the internal degrees of freedom of a chain are integrated out, many thousands of 
chains can be simulated within reasonable computer time. More recently, McCarthy, 
Lyubimov and Guenza introduced an approach based on the Ornstein-Zernike equations, 
in which a polymer is coarse-grained to a soft sphere on the scale of its radius of 
gyration [160J . At such a high level of coarse-graining typical maximum interaction 
energies are of the order of ksT, the thermal energy. 

All these approaches focus on reproducing the static structure of a polymer melt, 
with good to excellent agreement with the structure obtained from more detailed 
simulations. However, because such models ignore memory and/or entanglement effects, 
they do not have the characteristic slow dynamics of real polymer chains, and therefore 
cannot be used to predict the long-time rheology. 

One may be tempted to slow down the dynamics of such a super-coarse-grained 
system by using a large friction term. Lyubimov et al. recently introduced a method to 



estimate the friction on a polymer represented by a soft particle 161 . Their approach 



is based on an intermediate mapping of a chemically detailed chain on a freely-rotating 
chain model and using a hard-bead approximation to evaluate integrals involving the 
radial distribution function and dynamic structure factor. This yields a series of diffusion 
coefficients in good agreement with results of atomistic simulations of unentangled and 
slightly entangled polyethylene melts. No rheological predictions were made. 

The use of a single large friction term may be questioned, however, because in 
real entangled polymeric systems the frictions and random forces have memory of 
the configurations the system has gone through in the recent and sometimes even the 
distant past. A simple Brownian dynamics propagator with realistic mean forces and 
uncorrelated, fully random displacements without memory will not reproduce correct 
sequences of configurations of the retained coordinates. 

In recent years, our group has introduced and explored a new efficient framework to 



reintroduce memory of previous configurations in soft matter simulations 99 , 100 , 162 



166 . In this framework, called Responsive Particle Dynamics (RaPiD), configurations 
and forces are described by as few variables as possible. In case of polymer chains, 
similar to the works mentioned above, we usually restrict ourselves to three coordinates 
for the location of each centre of mass. To circumvent the complicated introduction 
of memory effects in friction forces and stochastic displacements (as in the Generalized 
Langevin Equation, Eq. (|2|) we introduce a relatively small set of additional dynamic 
variables, symbolically denoted as n{t), which keep track of the thermodynamic state 
of the eliminated coordinates for the given values of the retained coordinates. If the 
configuration of retained coordinates is suddenly changed (e.g. by a sudden change in 
the distance between the centres-of-mass of two polymers), the equilibrium values ng 
for the additional variables change too, but the relaxation to uq takes place over a finite 
time. This gives rise to strong transient forces in addition to the thermodynamic forces 
deriving from the potential of mean force. It is possible to show that, to lowest order, 
the transient forces may be thought of as originating from a penalty free energy function 
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which is quadratic in the deviations of the additional variables 163 



The particular additional variables and corresponding names that are used to 
describe the transient forces depend on the system under study. In some cases, such 
as for telechelic polymer networks, the variables are the number of bridging polymeric 
chains between nodes of the network 100, 165 . In other cases, more relevant to this 



review on polymer melts, the variables are the number of entanglements shared between 



a particular polymer chain and its immediate neighbours 162 



To some readers it will be rather surprising that a single particle model would 
produce rheological results that are most naturally described by assuming curvilinear 
tubes surrounding each chain which prohibit the chains to move in any way other then 
along the primitive paths of the tubes. A few qualitative remarks may be of help at 
this point 163 . First of all, the main action of the tubes is to conserve the prevailing 



configurations of the chains and their centres of mass in particular. This is exactly what 
the transient forces do as well. In a sense, the transient forces tie all centres of mass 
harmonically to their instantaneous positions. Secondly, beyond the Rouse time, i.e. as 
soon as the chains perform one-dimensional diffusive motions along their tube -axes, as 
a result of the Gaussian character of the tube-configurations the centres of mass of all 
chains perform simple random displacements. The only thing the tube does is to define 
the distribution of possible displacements and their averages in particular. Again this 
is what the transient forces do as well. 

The RaPiD method has demonstrated its usefulness by correctly reproducing the 
large-time linear and non-linear rheology of linear polymers 162 , telechelic polymers 



100 , 165 , star polymers 164 , and polymeric adhesives [166 . The next challenge is to 



link the few additional parameters introduced by the method to the chemical details of 
the system. In some cases, such as for the telechelic polymers 100,165 , such a link has 
already been successfully established. 



8. Conclusion 



The enormous range of time and length scales associated with the dynamics of entangled 
polymers precludes the use of one single computer simulation to calculate the large-scale 
rheological properties starting from the detailed chemical structure. Instead, some form 
of coarse-graining is necessary. In this review we have made a distinction between 
coarse-grained molecular dynamics (CGMD) and coarse-grained stochastic dynamics 
(CGSD). 

In CGMD a few atoms (up to five backbone carbon atoms) are lumped into 
relatively hard beads. Uncrossability of bonds is therefore automatically ensured. A 
large body of work has gone into deriving the effective interactions between the relatively 
hard beads, with iterative Boltzmann Inversion emerging as one of the most popular 
techniques. When coarse-graining the interactions, technically one should include 
additional friction and random forces, but fortunately in CGMD the friction is not 
dominant. If friction is ignored the dynamics will be too fast, usually by a factor 2 to 5. 
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Combined with the lower number of particles this leads to a significant computational 
speedup of order 100. Agreement with detailed atomistic simulations is often recovered 
by a suitable time scaling. Using these techniques it is possible to predict the rheology 
of melts of specific polymers up to one or two entanglement lengths. 

To reach a higher degree of coarse-graining, in CGSD many more atoms are lumped 
together (more than 10 backbone carbon atoms), leading to relatively soft beads. 
Friction and stochastic forces dominate the interactions. Unless counteractions are 
taken, the soft interactions cannot prevent crossing of the bonds between the coarse- 
grained particles. Possible counteractions include apphcation of the Twentanglement 
algorithm or prohibition of moves that lead to crossing of bonds. Using such techniques 
it is possible to predict the rheology of melts of specific polymers up to about 10 
entanglement lengths. 

To predict the rheology of polymers of longer polymers (say more than 10 
entanglement lengths), one is forced to make use of the tube model by obtaining 
entanglement characteristics through a primitive path analysis and simulating a 
primitive chain network. Primitive path analysis algorithms allow us to determine the 
tube diameter and entanglement molecular weight from atomistically detailed or CGMD 
simulations. The entanglement molecular weight so obtained is often in good agreement 
with the entanglement molecular weight derived from the experimental platcaii modulus. 
However, predictions of timc-dcpcndcnt properties rely on details of the interpretation of 
reptation theory, and often seem to be determined by a somewhat different entanglement 
molecular weight. A prerequisite for future computational developments is to have 
available a better understanding of the relation between the two apparently different 
entanglement molecular weights and their relations to the effects of chain stiffness, finite 
size fiuctuations and tube renewal mechanisms. This information needs to be reliably 
predicted from detailed, chemically realistic models. Improvements in understanding 
will be used to develop more realistic slip-link models. In this respect it is important 
to note that already now slip-link models are very effective and fast in deriving the 
rheological properties of polymer melts. 

Finally we discussed the use of super-coarse-grained models, in which a full polymer 
is represented by a single particle, which may be used to predict fiow and rheology on 
industrially relevant large time and length scales. Because coarse-graining has been 
taken to its limit, conservative interactions between full polymer chains are extremely 
soft and memory or entanglement effects are prominent. To efficiently introduce memory 
effects we have developed the Responsive Particle Dynamics method. 

In summary, our search for the link between chemical structure and large time 
dynamics has led to many innovations during the last two decades. Computer 
simulations will continue to smoothen the bumpy road from chemistry to rheology. 
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